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I .  INTRODUCTION 


In  this  report  optimal  designs  for  weighted  least  squares  and 
maximum  likelihood  estimation  of  the  two-parameter  logistic  function 
are  constructed^  First,  some  properties  of  the  logistic  function  are 
discussed  along  with  techniques  for  its  estimation.  Then,  four 
criteria  for  optimality  are  defined  and  the  corresponding  optimal 
designs  are  constructed.  Finally,  some  practical  considerations  in 
the  implementation  of  these  designs  are  mentioned. 


1 

A  major  portion  of  the  research  described  in  this  report  is  based 
on  the  senior  author's  Masters  paper  [3]  at  the  Pennsylvania  State 
University.  That  paper,  written  under  the  direction  of  Dr.  James  L. 
Rosenberger  of  the  Statistics  Department,  has  been  issued  by  that  depart¬ 
ment  as  Technical  Report  No.  33  (August  1978). 
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II.  THE  LOGISTIC  FUNCTION 


The  logistic  function  has  the  form 

P  “  f(x,  J3)  “  [1  +  exp(-x'ji)]  1  (1) 

where  x  *  (1,  Xj_,  X2>...,Xjc)'  and  =  (3Q,  32»  *  *  *  ’V  *  P  takes  on 

values  in  the  interval  (0,  1),  8  is  the  vector  of  parameters  and  x  is 
the  vector  of  independent  (predictor)  variables.  The  first  element  in 
x,  a  "dummy  variable",  is  included  to  provide  for  the  estimation  of  an 
intercept. 

This  function  is  often  used  to  describe  the  relationship  between 
the  vector  x  and  the  probability  of  a  certain  response,  where  the  response 
is  dichotomous.  For  example,  the  U.S.  Navy's  impact  acceleration  research 
program  being  conducted  by  the  Naval  Aerospace  Medical  Research  Laboratory 
(NAMRL)  Detachment  is  concerned  with  the  relationships  between  various 
dynamic  quantities  (e.g.,  peak  head  linear  acceleration  and  peak  head 
angular  acceleration)  and  the  probability  of  injury.  Here,  x  represents 
the  values  of  the  dynamic  measurements  and  P  is  the  corresponding  probability 
of  injury. 

A.  ESTIMATION  OF  THE  LOGISTIC  FUNCTION 

Letting  the  n  observations  in  an  experiment  be  distinguished  with 
the  subscript  i  =  1,  2,...,n,  the  observed  probability,  denoted  p^,  is  given 
by 

Pi  “  Pi  +  £i  =  f^-i’  +  ei 
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where  denotes  the  error  term.  Because  of  the  binomial  response  (e.g., 
injury  or  noninjury)  pi  is  either  1  or  0.  Therefore,  is  either  Q±  -  1  -  P 
or  with  corresponding  probabilities  P^  and  respectively.  From  this  it 
follows  that 
E(et)  -  0 
Var  (e±)  *  PiQi 

Since  the  error  variances  are  not  necessarily  equal  for  different  observations, 
it  is  appropriate  to  use  weighted  least  squares  to  estimate 

It  can  be  shown  that  the  weighted  least  squares  estimate  of  j5  is  given 
by 

b  -  (X'HX)"1  X'Hy  (2) 

where  b  denotes  the  estimated  parameter  vector 

J>  m  (Bq,  6^,  • • • f 3^) ’i  (3) 

X  denotes  the  design  matrix 


H  denotes  the  diagonal  weight  matrix 


and  £  denotes  the  vector  of  "working  observations" 


z  =  (yx»  y2»*”»yn)'  (6) 

where  y  =  _1 _  [p  -  P  +  P  Q  ln(P  /Q  )]. 

x  po  x  x  ii  ii 
iwi 

For  a  general  discussion  of  weighted  least  squares,  including  a  derivation 

of  equation  (2),  see  Draper  and  Smith  [1].  Equation  (2)  is  identical  to 

that  which  would  be  obtained  by  maximum  likelihood  estimation.  For  a  more 

complete  mathematical  development  of  the  estimation  process  see  Walker 

1 

and  Duncan  [8]  or  Smith  [7]. 

Note  that  since  H  and  ^  are  functions  of  the  probabilities 

Pi  *  [1  +  exp(-x'jS)]  1,  the  estimate  t>  in  equation  (2)  must  be  obtained 

A 

iteratively  by  using  the  values  P  =  [1  +  exp(-x'b)  ]  from  the  previous 
iteration.  Because  each  observation  (of  injury  or  noninjury)  corresponds 


1 

Equation  (2)  and  definitions  (3)  through  (6)  are  mathematically 
equivalent  to  the  corresponding  formulas  given  in  the  references  [7,  8]. 
However,  changes  in  notation  were  made  so  that  formulas  could  be  presented 
as  they  might  be  found  in  a  standard  textbook  discussion  of  weighted  least 
squares  (such  as  Draper  and  Smith  [1].) 


to  a  value  of  =  1  or  =  0»  the  actual  observations  cannot 

be  used  as  Initial  estimates.  This  Is  not  a  major  difficulty,  and  may 
be  overcome  by  using  Initial  estimates  obtained  by  fitting  a  discriminant 
function.  (See  Jones  [2].)  An  alternate  recursive  procedure  (Walker  and 
Duncan  [8])  requires  that  initial  estimates  of  the  parameter  vector  and 
its  covariance  matrix  be  specified. 

The  asymptotic  covariance  matrix  of  b  is 

Var(b)  =  (X'HX)'1  .  (7) 

In  practical  applications,  the  covariance  matrix  must  be  estimated  by 
substituting 

A 

H  = 


for  H  into  (7) . 

B.  THE  TWO-PARAMETER  CASE 

Since  this  report  discusses  optimal  designs  for  the  two-parameter 
case  (k  ■  1),  it  will  be  convenient  to  introduce  some  additional  notation 
and  terminology  for  that  special  case.  When  k  =  1,  equation  (1)  is  a 
sigmoid  curve  of  the  form 

P  -  {1  +  exp[-(0o  +  ^x)]}"1  (8) 


P1Q1 


P  Q 
n  n 
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and  for  >  0,  P  has  asymptotic  limits  0  and  1  as  x  approaches  -°°  and  <*> 
respectively.  (See  Figure  1.) 

As  an  outgrowth  of  the  early  biological  applications  of  the  logistic 
function,  x  is  often  referred  to  as  the  "dose"  and  the  level  of  x  correspond¬ 
ing  to  a  probability  P  of  response  is  denoted  where  LD  refers  to 

"lethal  dose."  For  example,  a  dose  of  amount  LD75  would  result  in  a 
response  with  probability  .75.  Designs  for  estimating  the  two-parameter 
logistic  function  are  characterized  in  terms  of  LD  levels  at  which  obser¬ 
vations  are  taken.  For  example,  a  typical  design  may  allocate  one  quarter 
of  the  experimental  units  to  LD^q,  one  half  of  the  experimental  units  to 
LD^q  and  one  quarter  of  the  experimental  units  to  LD^q. 
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III.  OPTIMAL  DESIGNS 


In  general,  an  optimal  design  is  one  which  "minimizes"  the  covariance 
matrix  (X'HX)  The  meaning  of  the  word  minimize,  when  applied  to  a 
matrix,  is  not  obvious;  a  number  of  functionals  of  a  covariance  matrix 
have  been  proposed  as  criteria  for  minimization.  For  each  criterion  there 
is  a  corresponding  optimal  design.  (See  Kiefer  [4].)  These  include: 

1.  D-optimality .  Minimize  by  choice  of  design,  the  determinant 
of  (X'HX)”1,  denoted  | (X'HX)-1 | . 

2.  A-optimality .  Minimize  by  choice  of  design,  the  trace  of 
(X'HX)-1,  denoted  tr ( (X'HX)-1) . 

3.  E-optimality .  Minimize,  by  choice  of  design,  the  maximum 
eigenvalue  of  (X'HX)  1. 

4.  G-optimality .  Minimize,  by  choice  of  design,  the  maximum 
variance  (over  all  dose  levels  x^)  of  the  predicted 
var(P^),  where 

PQ  =  U  +  exp [- (3q  +  BjXq)]}-1  . 

These  criteria  are  now  examined  more  closely.  Clearly,  all  of  the 
information  in  an  r  x  r  covariance  matrix  needed  by  the  D,  A  and  E-optimality 
criteria  is  available  in  its  eigenvalues,  denoted  A-p  X2,...,Xr,  since 

.  _1  r  _ l  r 

|  (X'HX)  |  =  II  A  and  tr[  (X'HX)  ]  =  I  A..  These  three  criteria  can 
i=l  i-1 

therefore  be  expressed  (for  the  two-parameter  case): 
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1. 

D-optimality . 

Minimize 

<X1 

*  *2). 

2. 

A-optimality . 

Minimize 

<A1 

+  A2) . 

3. 

E-optimality. 

Minimize 

max 

^ 2^  * 

The  determinant  of  (X'HX)  ,  I  (X'HX)  =  A^  •  a  general 

measure  of  the  precision  of  the  estimates  3q  and  3^.  Note  that  | (X’HX)  | 

A  A  A  A  n 

=  var (3g)var (3^)  -  [(cov(Bq,  3^)]  .  It  will  be  shown  that  the  D-optimal 
design  is  invariant  to  the  choice  of  parameters  or  units  of  measurement  of 
the  dose. 

The  trace  of  (X'HX)  \  tr [ (X'HX)-’'-]  =  A^  +  A^,  is  not  as  general  a 

i  — 1 1 

measure  of  precision  as  | (X'HX)  | ,  since  the  trace  criterion  ignores 

A  A  _  ^ 

information  about  the  covariance  of  Sq  and  3^.  Note  that  tr [ (X'HX)  ]  = 

A  A 

var(3Q)  +  var(B^).  Furthermore,  the  A-optimal  design  is  not  invariant 
to  the  choice  of  parameters  or  units  of  measurement  of  the  dose.  The 
practical  implications  of  this  fact  will  be  discussed  in  Section  IV. 

The  maximum  eigenvalue  of  (X'HX)  \  max  (A^,  A^),  is  the  variance 

A.  /\ 

of  that  linear  combination  of  Bq  and  3-^  which  is  least  precisely  estimated 
by  the  design.  Thus,  E-optimality  is  a  minimax  criterion  on  the  variance 
of  all  linear  combinations  of  the  estimated  parameters.  Since  the  set  of 
eigenvalues  gives  all  of  the  information  contained  in  (X'HX)  ,  some  infor¬ 
mation  is  lost  in  ignoring  the  smaller  eigenvalue.  Like  the  A-optimal 
design,  the  E-optimal  design  is  not  invariant  to  the  choice  of  parameters 
or  units  of  measurement  of  the  dose. 

In  considering  the  G-optimal  design,  it  should  be  noted  that  the 
variance  of  a  predicted  value  is  a  function  not  only  of  (X.'HX)  \  but  also 


of  Xq,  the  dose  at  which  the  probability  is  being  predicted.  Furthermore, 


the  G-optimality  criterion  is  a  mioimax  criterion  on  the  variance  of  a 
non-linear  function  of  8^  and  8^,  namely  the  logistic  function  (8).  Like 
the  D-optimal  design,  the  G-optimal  design  is  invariant  to  the  choice  of 
parameters  or  units  of  measurement  of  the  dose.  It  can  be  shown  (Kiefer 
and  Wolfowitz  [5])  that  the  D-optimality  and  G-optimality  criteria  are 
equivalent  for  ordinary  least  squares  estimation.  (Recall  that  weighted 
least  squares  estimation  is  preferred  for  the  logistic  function.) 

The  four  criteria  discussed  above  are  applicable  regardless  of  the 
number  of  levels  of  x  at  which  observations  are  taken.  However,  only  two- 
point  designs  will  be  considered  in  this  report  for  the  following  reasons. 
Sibson  and  Kenney  {6}  have  shown  that  for  ordinary  least  squares  estimation, 
the  I>-optimal  (and  equivalently  the  G-optimal)  designs  for  estimating  r-order 
polynomial  functions  are  (r  +  l)-point  designs.  The  logistic  function  (8) 
can  be  written  as  a  first  order  polynomial  as  follows: 

ln[P/(l  -  P)]  =  Bq  +  BjX  .  (9) 

The  transformation  by  which  the  function  (8)  was  linearized,  g(P)  = 
ln[P/(l  -  P)],  is  often  called  the  "logit"  or  "log  odds  ratio." 

Although  weighted  least  squares  estimation  is  being  used,  it 
seems  reasonable  to  believe  that  the  optimal  designs  are  two-point  designs. 
Furthermore,  since  the  logistic  function  is  symmetric  with  respect  to  two 
rotations  about  x  =  LD,-q  and  P  =  .50  [that  is,  f(LD^QQp)  =  1  -  f  (LD^qqq_P)  )  ] 
only  symmetric  designs  about  LD^q  with  equal  allocation  to  both  design  points 
will  be  considered. 

Letting  Q'  ■  1  -  P',  the  design  points  will  be  denoted  x^  =  LD^QQpi 
and  X2  =  w^ere  P'  y  -50.  Thus,  P1  and  Q'  are  the  values  of  the 
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logistic  function  evaluated  at  design  points  x^  and  x2  respectively. 

The  sample  size,  denoted  n,  is  assumed  to  be  sufficiently  large  to  use 
asymptotic  theory. 

Hence  the  symmetric  two-point  designs  can  be  indexed  simply  by  P'. 
If,  for  example,  P'  =  .75,  the  design  of  the  experiment  would  be  to 
randomly  allocate  one  half  of  the  observations  at  x^  ■  LD75  an<*  one  half 
of  the  observations  at  =  LD25.  (See  Figure  2.) 

Using  the  aforementioned  notation  and  restrictions,  X'HX  can  be 
reduced  to 


2 

X'HX  =  (n/2)P'Q'  i 

-  i  X1  +  x2 

A.  CONSTRUCTION  OF  D-OPTIMAL  DESIGNS 

The  problem  is  to  minimize  | (X'HX)  over  choices  of  design.  This 

is  equivalent  to  maximizing  |X’HX | ,  since  for  any  nonsingular  matrix  A, 

I A  1 1  =  1/ 1 A | .  Thus,  the  D-optimality  criterion  can  be  written  max|X'Hx|. 

P' 

Now  from  (10), 

|X'HX|  =  [(n/2)P'Q'(x1  -  x,,)]2.  (11) 

But  from  (8)  or  (9), 

x!  =  [ ln(P ' /Q * )  -  Bq]/&i  and  x^  =  [ln(Q'/P')  -  (12) 

Substituting  these  equalities  into  (11), 


Xj^  +  X2 


(10) 
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|X'HX|  =  [nP'Q'liUP’/Q')/^]2. 


(13) 


So  max  | X ' HX |  *  maxtnP'Q'lnCP'/Q'J/B-,  ]2. 
p.  p.  1 

To  maximize,  set  the  first  partial  derivative  with  respect  to  P'  equal 
to  zero.  The  solution  to  this  can  be  numerically  obtained  as  P*  =  .824, 

Q'  =  .176. 

Thus,  the  two-point,  symmetric  about  LDjq,  D-optimal  design  is  to 
allocate  half  of  the  experimental  units  to  LDq2  4  and  half  of  the  experi¬ 
mental  units  to  LD, ,  ,.  Note  that  the  optimal  design  is  invariant  to  the 
1  / .  o 

choice  of  parameters,  3q  and  Bi»  or  units  of  measurement  of  the  dose,  x. 

B.  CONSTRUCTION  OF  A-OPTIMAL  DESIGNS 

-1 

The  trace  criterion  specifies  that  tr [ (X'HX)  ]  be  minimized  over 
choices  of  design.  From  (10)  and  (12)  it  can  be  shown  (see  Kalish  [3]) 
that 

tr [  (X'ffiC)-1]  =  1/nP'Q'  +  (3g  +  eJ)/{nPQ[ln(P7Q')]}2 

To  minimize,  set  the  first  partial  derivative  with  respect  to  P'  equal  to 

zero.  The  solution  is  not  invariant  to  the  choice  of  8q  and  In  fact, 

the  A-optimal  design  depends  on  the  sum  of  the  squares  of  the  parameters 
2  2 

(i.e.,  on  8g  +  3^).  A-optimal  designs  were  constructed  for  several  choices 

of  parameters,  using  a  computerized  numerical  search  approach  (see  Kalish  [3]). 

Figure  3  shows  the  relationship  between  the  A-optimal  design  (indexed  by  P') 

2  2  2  2 

and  8q  +  3^.  Note  that  as  8q  +  3^  gets  larger,  the  design  points  move  out- 

2  2 

ward  from  1°  fact,  their  asymptotic  limits  (as  Bq  +  8^  approaches  °°) 

are  ^  and  LDqq  3* 
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C.  CONSTRUCTION  OF  E-OPTIMAL  DESIGNS 


For  this  criterion,  the  maximum  eigenvalue  of  (X'HX)  is  minimized 
over  choices  of  design.  Recall  that  the  eigenvalues  of  a  k  x  k  matrix, 

M,  are  the  k  roots  of  the  equation  |m  -  Al^ |  =  0,  where  1^  denotes  the 
identity  matrix  of  order  k.  It  can  be  shown  (see  Kalish  [3])  that  the 
maximum  eigenvalue  is  given  by 

00  +  f*i  +  Un(P'/Q'))2  +  /TbJ  +  Si  +  (ln(P’  /Q'  ))^T^-[23^1n(P’  /Q’ )] 2 
maxCAj^,  A2)  =  _ 

2  nP ' Q ' (In (P'/Q'))2 

The  problem  of  minimizing  this  equation  with  respect  to  P'  was  approached 
analytically  but  the  calculations  became  intractable.  Thus,  a  computerized 
numerical  search  was  used  (see  Kalish  [3]).  As  with  the  A-optimal  designs, 
the  E-optimal  design  points  move  outward  towards  LD^  ^  and  LDqq  ^  as  the 
parameters  get  larger.  In  this  case,  however,  P'  is  not  a  monotone  function 
of  0q  +  ej. 

D.  CONSTRUCTION  OF  G-OPTIMAL  DESIGNS 

Recall  that  G-optimality  is  a  minimax  criterion  on  the  variance  of  the 
predicted  values.  As  before,  denote  the  probabilities  corresponding  to 
design  points  x^  «  LD100p,  and  x2  =  LDioO(l-P')  as  P'  and  Q'  resPectively • 

A 

Furthermore,  refer  to  as  the  predicted  value  corresponding  to  x^.  Using 
this  notation,  the  G-optimality  criterion  can  be  expressed  min  max[var (P.) ] . 

pi  S  0 

*  S  p0 

Note  that  var(PQ)  cannot  be  calculated  directly  since  Pq  is  a  nonlinear 

function  of  8q  and  8^.  However,  using  a  Taylor  series  expansion,  it  can  be 

A 

shown  that  the  asymptotic  variance  of  Pq  is 
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var(PQ) 


(*A  “  Xl)  ^  (^0  “  *o)  T  t _  O  O 

01  02  [e'L0/(l  +  e-L0)2]2, 


(n/2)P'Q'(Xl  -  x2)‘ 


where  Lq  *  1ii[Pq/(1  -  Pq)]  “  Bq  +  B^Xq.  This  can  be  further  simplified 
to  be  written  only  in  terms  of  n,  Pq,  P'  (and  Qq  =  1  -  Pq,  Q'  -  1  -  P')  as 
a  [ln(PoQVQoP')]2  +  [^(PqP'/QoQ')]2  2 


Var (Pq) 


(n/2)P’Q'[ln(P2/Q2)]2 


The  fact  that  var(PQ>  can  be  expressed  as  a  function  only  of  Pq  and  the 
design  (indexed  by  n  and  P)  shows  that  the  G-optimal  design  must  be 
independent  of  the  parameters  8q  and  6^.  Thus,  the  G-optimality  criterion 
is  invariant  to  the  choice  of  parameters  or  units  of  measurement  of  the 
dose. 

An  analytic  approach  to  the  minimax  problem  was  attempted  but  again 
the  calculations  became  intractable.  Thus,  a  computerized  numerical  search 
was  used  (see  Kalish  [3]).  It  was  found  that  the  G-optimal  design  has 
design  points  at  x^  ■  LD^g  g  and  x2  =  LD  23  2" 


IV.  PRACTICAL  CONSIDERATIONS 


The  practical  applications  of  the  optimal  designs  discussed  here 
require  several  considerations.  One  problem  lies  in  the  fact  that  the 
specification  of  an  optimal  design  is  done  in  terms  of  LD  levels,  yet  in 
the  design  stage  of  an  experiment,  the  LD  levels  are  typically  not  known. 
Obviously,  if  a  substantial  amount  of  prior  knowledge  is  available  before 
experimentation,  this  can  be  used  to  give  estimates  of  LD  levels  and  of 
parameter  values  in  order  to  approximate  an  optimal  design. 

If  little  or  no  prior  information  is  available,  a  small  pre-study 
experiment  might  be  conducted  wherein  experimental  units  are  allocated 
to  some  "reasonable"  range  of  doses  expected  to  cover  the  LD^,.  to  LD^ 
levels.  Data  from  the  pre-study  can  then  be  used  to  estimate  parameters 
and  optimal  design  points.  Of  course,  the  pre-study  data  can  later 
augment  the  primary  experimental  data  to  obtain  final  estimates. 

Another  practical  problem  to  be  resolved  is  which  optimality 
criterion  to  use.  It  has  already  been  mentioned  that  the  D-optimal  and 
G-optimal  designs  are  invariant  to  choice  of  parameters,  while  the  A- 
optimal  and  E-optimal  designs  are  not.  Since  the  values  of  a  parameter 
are  interpreted  in  the  same  units  as  x,  any  change  in  the  scale  of  measure¬ 
ment  of  x  would  result  in  different  A-optimal  or  E-optimal  designs.  For 
example,  if  velocity  (x)  were  measured  in  meters  per  second,  the  A  and 
E-optimal  designs  would  differ  from  the  A  and  E-optimal  designs,  respectively, 
if  velocity  were  measured  in  centimeters  per  second. 

Furthermore,  it  has  been  noted  that  the  A  and  E-optimality  criteria 
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do  not  use  all  of  the  information  contained  in  (X'HX)  The  trace 

A  A 

criterion  considers  only  the  diagonal  elements,  varCB^)  and  var(S^),  of 

A  A 

the  covariance  matrix  but  not  the  off-diagonal  elements,  cov(3q»  3^)- 
The  maximum  eigenvalue  criterion  utilizes  the  variance  of  only  one  of 

A  A 

two  orthogonal  linear  combinations  of  3q  and  B^.  For  these  reasons, 
the  D-optimality  and  G-optimality  criteria  seem  superior  to  the  A- 
optimality  and  E-optimality  criteria. 

Assuming,  then,  that  an  experimenter  has  yet  to  choose  between  D 
and  G-optimality,  the  following  observation  makes  the  decision  seem  less 
critical.  Consider  plots  of  these  two  criteria  versus  P'.  Graphs  of 
| (X'HX)  |  versus  P'  and  max  var(PQ)  versus  P'  are  displayed  in  Figures  4 
and  5.  It  can  be  seen  that  both  curves  are  fairly  flat  for  values  of  P' 
between  .75  and  .85.  Due  to  the  flatness  of  these  curves  in  the  region 
of  the  optimal  designs  (i.e.,  P'  =  .824  for  D-optimality  and  P'  *  .768 
for  G-optimality),  one  can  "miss"  the  optimal  designs  and  still  not 
sacrifice  much  efficiency  in  the  estimation  of  the  function.  In  fact, 
one  can  "almost"  achieve  optimality  for  both  criteria  simultaneously. 


V.  SUMMARY 


In  this  report  optimal  designs  were  constructed  for  estimation  of 
the  two-parameter  logistic  function.  In  particular,  four  criteria  for 
optimality  were  used:  D,  A,  E  and  G-optimality.  It  was  shown  that  the 
D  and  G-optimality  criteria  are  invariant  to  changes  in  scale  or  units 
of  measurement  of  the  independent  variable.  In  addition,  the  A  and  E- 
optimality  criteria  ignore  some  of  the  information  available  in  the 
covariance  matrix.  For  these  reasons,  the  D  and  G-optimality  criteria 
seem  superior  to  the  A  and  E-optimality  criteria. 

Since  no  criterion  can  be  applied  exactly  in  a  real  life  setting, 
the  problem  of  approximating  an  optimal  design  was  briefly  discussed. 
Two  future  technical  reports  will  discuss  an  extension  of  this  idea: 
that  of  augmenting  an  existing  design  one  point  at  a  time  using  an 
optimality  criteria  closely  related  to  D-optimality . 
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